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Abstract 

Primary GCR interact with the Earth's atmosphere originating atmospheric 
showers, thus giving rise to fluxes of secondary particles in the atmosphere. 
Electromagnetic and hadronic interactions interplay in the production of these 
particles, whose detection is performed by means of complementary techniques 
in different energy ranges and at different depths in the atmosphere, down to the 
Earth's surface. 

Monte Carlo codes are essential calculation tools which can describe the 
complexity of the physics of these phenomena, thus allowing the analysis of 
experimental data. However, these codes are affected by important uncertainties, 
concerning, in particular, hadronic physics at high energy. In this paper we shall 
report some results concerning inclusive particle fluxes and atmospheric shower 
properties as obtained using the FLUKA transport and interaction code. Some 
emphasis will also be given to the validation of the physics models of FLUKA 
involved in these calculations. 
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1 Introduction: EAS development and detection 

GCR interactions in the Earth's atmosphere, induced by particles with energies high 
enough (E > 10^'' — 10^^ eV), may originate Extended Air Showers (EAS). Both the EM 
and the component of EAS can be detected and used to infer properties concerning 
GCR primary spectrum energy and composition. At present, a special emphasis is put 
in the investigation of Very High and Ultra High Energy Cosmic Rays (VHECR and 
UHECR), through dedicated experiments devoted to the detection of EAS originated by 
these particles, due to the still open questions concerning the shape of the spectrum (2"°' 
knee, ankle), the composition (7, p, heavy ions), and the origin of cosmic rays (galactic 
/ extragalactic), their transport, acceleration and the GZK cut-off. 
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Information on GCR primary flux and composition are inferred from tlie experimental 
observables by using Monte Carlo (MC) simulations. In particular, while the fluorescence 
technique, used by many experiments, allows with some uncertainties to estimate the EM 
energy deposited in the atmosphere, the total primary energies can be obtained only after 
estimating the missing energy carried by other shower components (/i, nucleons...), and 
this requires indeed MC simulations. Even the search for particular estimators, such as 
the S variable, allows to obtain primary energy information in a way nearly independent 
of composition [1], but needs MC for calibrations. 

In particular, in the AUGER experiment |2] information on primary energies are in- 
ferred by hybrid measurements and the fluorescence technique as far as the EM sector 
is concerned (with an uncertainty amounting to ~ 20%), whereas the number of /xs at 
the Earth's surface, N^, is still related to total primary energies by means of MC simu- 
lations. On the other hand, in the KASCADE and KAS CADE- Grande experiments [3] 
information on primary energies are inferred from measures of fi and e with the aid of MC 
shower simulations. Both in AUGER and in KASCADE / KASCADE-Grande, as well 
as in many other EAS experiments, MCs are heavily needed to analyse data in terms of 
GCR mass composition. At present there are however still important uncertainties, since 
different MC models produce important differences on the interpretation of the same 
data. The critical ingredients in these MC codes are the physical models for hadronic 
and nuclear interaction, and their implementation. 

In this work we summarize some of the results which can be obtained by MC simula- 
tion of EAS, giving examples taken from the use of the FLUKA code. 

2 MC simulation of EAS 

EAS simulations can be performed by means of MC codes, which include hadronic and 
EM modules. This method leads, in a straightforward way, to account for fluctuations in 
the evolution of different showers initiated by primaries of the same type and energy. 

The GCR primary spectrum covers a wide energy range, extending on many orders 
of magnitude. This fact implies that it is very difficult to have in a MC code a single 
hadronic model able to describe primary and secondary hadron interactions in the Earth's 
atmosphere. Therefore in many cases different models, differing according to the energy 
range under study, have to be merged together. 

At present the most diffused MC package for GCR induced shower simulation is 
CORSIKA [1]. As far as its hadronic sector is concerned, it distinguishes between high 
energy (E > 100 - 200 GeV) and low energy (E < 100 - 200 GeV) hadronic models. 
In the present version, CORSIKA offers the choice, at high energy, among models like 
QGSJET, QGSJET-II, SIBYLL, DPMJET-2.55, NEXUS, EPOS, whereas, at low energy, 
GHEISHA, UrQMD, FLUKA are available. Other packages for GCR EAS simulation also 
exist. We just mention, among the others, the AIRES [5j and the Cosmos [U] codes. 

In this paper, we focus instead on the performances of the FLUKA multipurpose 
code [7j. FLUKA is a fully integrated, high precision, transport and interaction code, 
without the need of invoking other external packages. Besides its use in the CORSIKA 
package at energies < 100 - 200 GeV, FLUKA has recently been used also as a standalone 
code, complemented by the DPMJET-III code [8], for the simulation of EAS at energies 
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< 10^^ eV. FLUKA simulations can, in principle, lead to results different from the ones 
made by FLUKA+CORSIKA, due to the implementation of different EM and hadronic 
models. FLUKA contains a highly accurate EM interaction model [9l , not included in 
CORSIKA. As far as the hadronic sector is concerned, the following models are available 
in FLUKA: 

- low energy (E < 20 TeV) h + A reactions can be simulated by means of the PEANUT 
(PreEquilibrium Approach to NUclear Thermalization) module, originally worked out 
for energies < tens of GeV [H], and recently improved and extended to higher ener- 
gies |12] . This extension has superseded a previously existing DPM + INC model, working 
at E > 5 GeV. 

- high energy (E > 20 TeV) h + A reactions can be described by means of the interface 
to the DPMJET-2.5 or III codes. 

- A + A reactions can be simulated by means of an interface to the RQMD2.4 code, 
written in Frankfurt and modified for insertion in FLUKA [13], for E < 5 GeV/ A, and by 
means of the interface to the DPMJET-2.5 or III codes for E > 5 GeV/A. Alternatively, 
the superposition model is also available to roughly describe heavy-ion interactions. 

Constraints on theoretical models for particle interactions can (or will) be obtained 
by data collected at collider experiments, but, unfortunately, diffractive cross-sections are 
difficult to be measured at colliders. Furthermore, present and near future accelerators 
do not and will not allow to test the physics of particles at the highest GCR energies 
(E > 1.5 10^^ eV, GCR spectrum tail). Thus, MC models used in the analysis of EAS 
data are still affected by important uncertainties, concerning in particular the hadronic 
sector. Other constraints/checks on low-energy models, besides the ones from accelerator 
data, can come from astrophysical measurements of inclusive particle fluxes at different 
atmospheric depths, as performed by balloon-borne experiments and satellites. We em- 
phasize that for an accurate description of VHECR and UHECR EAS also low-energy 
physics is important, e.g. the larger is the distance from the shower core and the shower 
inclination, the stronger is the required accuracy of the adopted low energy models. 

2.1 Validation of the FLUKA models at low energies 

Validation of the FLUKA models at low energies has been performed by means of com- 
parisons of their theoretical predictions with data collected at the accelerators. Data 
concerning inclusive particle fluxes in the atmosphere down to the Earth's surface, have 
been then used to assess the FLUKA performance on cosmic ray physics. 

While the first issue can be addressed by using the version of the code officially 
distributed and available on the web, the second one requires the inclusion of further (ge- 
ometrical, geomagnetical, astrophysical) assumptions and information. These elements 
have been inserted by means of the FLUKA user routines. 

Validation of FLUKA using accelerator data Data recently collected by the HARP 
and the NA49 experiments have been nicely reproduced by means of FLUKA [10] , Many 
others validations of the code, referred to previously published data, can be found in 
literature. Furthermore, data from the SPS collider, Tevatron and RHIC have been used 
to successfully validate the DPMJET2.5 and III codes [T4] . 
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Figure 1: Calculated inclusive fluxes of /i~ [left panel) and /x"*" {right panel) at the top of 
Mt. Norikura (~ 2700 m a.s.l.) as obtained by FLUKA simulations (open symbols) vs. 
exp. data [l5] from BESS'99 (full symbols). Other examples of FLUKA benchmarks in 
GCR physics can be found in [16] . 



Figure 2: Calculated inclusive fluxes of fi^ {left panels) and //^ {right panels) at atmo- 
spheric depths of 11.9 g/cm^ {odd panels) and 26.4 g/cm^ {even panels), respectively, as 
obtained by FLUKA simulations (open symbols) vs. exp. data [17] from BESS 2001 (full 
symbols) . 



Validation of FLUKA using astrophysical data from satellites and balloon- 
borne experiments Inclusive particle fluxes in the Earth's upper and lower atmo- 
sphere, induced by GCR propagation and interactions, have been measured by means of 
satellites and balloon-borne experiments such as BESS. Different sets of data, taken in 
different atmospheric and geomagnetic conditions, have been analyzed and reproduced 
by FLUKA, under proper assumptions concerning GCR primary spectrum, solar modula- 
tion, atmospheric and geomagnetic models. As examples, data concerning /i collected at 
the top (~ 2700 m a.s.l.) of Mt. Norikura (geomagnetic cut-off ~ 11.2 GV) in 1999, and 
at two different atmospheric depths (11.9 and 26.4 g/cm^) in a balloon flight campaign 
over Ft. Sumner (geomagnetic cut-off ~ 4.2 GV) in 2001, are presented in fig. [T] and 
[21 respectively, and compared to the theoretical fluxes expected on the basis of FLUKA 
simulations. The agreement between the experimental and theoretical results is quite 
nice and increases at increasing energies. Furthermore, atmospheric /x charge ratios in 
the energy range from ~ 20 GeV up to ~ 3 TeV have been detected by the L3+C expe- 
riment at CERN (~ 450 m a.s.l.), that has provided also /i arrival directions. FLUKA 
has successfully reproduced also these data [TH] . 

Proton and heavier ion inclusive fluxes have also been obtained by FLUKA theoretical 
simulations and compared to experimental data available from BESS. In these cases, the 
results of the simulations turn out to be very sensitive to a correct treatment of the 
geomagnetic cut-off, of the transport of particles in the upper layers of the atmosphere 
(trapping of particles and their recirculation), and to the accuracy of the theoretical 
modelling of low-energy A-A interactions. 

2.2 Predictions of the FLUKA models at high energies: EAS 

In the following we shall report about some tests of EAS simulation at energies 10^^ - 
10^^ eV, performed by means of FLUKA, complemented by DPMJET-llI for the descrip- 
tion of A-A reactions at E > 5 GeV/A and h-A reactions at E > 20 TeV. Vertical showers, 
induced by primaries of different mass {p and Fe ions), have been considered. 

The Vmax position, i.e. the atmospheric depth where the shower reaches its maximum 
development (in terms of EM energy deposition), is a good probe of GCR composition. 
In fact, since the inelastic scattering cross-section for ;3-Air interactions is lower than the 
A-Air ones, EAS induced by Fe primaries have, on the average, vertical profile maxima 
shifted towards lower atmospheric depths than p induced showers, e and fluences as a 
function of the atmospheric depth p for p and Fe induced showers at 10^^ and 10^^ eV, 
as calculated with FLUKA + DPMJET-lIl, are shown in ref. [TU]. In general, different 
codes, based on different models, due to different modelings of hadronic and EM processes, 
give different results for Amax- Anyway, to disentangle the effects of different models and 
different primary composition, also more specific variables can be considered. Among 
the others, the study of Amax fluctuations [E] for primaries at a fixed energy has been 
proposed for this purpose. As an example, Amax fluctuations for Fe induced showers have 
been computed from our simulations, considering the e, p, charged hadron and neutral 
hadron vertical shower profiles, and are shown in the central and right panels of fig. [3l In 
general, the average muon Xmax is located at larger ps with respect to the electron one, 
since /is penetrate more deeply. The charged hadrons undergo strong absorption effects, 
thus their Xmax are located at lower ps with respect to the ones of neutral hadrons. 
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Figure 3: Left panel: fluctuations in the distribution of e X^ax for a few hundred vertical 
showers induced by p {dot-dashed green line) and Fe ions [dotted red line) at 10^^ eV 
energy. Central and Right panels: fluctuations in the distribution of e X^ax [central 
upper panel), /i X^ax [central lower panel) , charged hadron X^ax [right upper panel) and 
neutral hadron Xmax [right lower panel) for a few hundred vertical showers induced by 
Fe ions at 10^^ eV [dotted red line) and 10^^ eV [dashed blue line) energies. 



Increasing energies lead, in general, to a shift of X^ax to higher ps for all considered 
profiles, even if, in all cases, the fluctuation profiles at 10^^ and 10^^ eV overlap. That 
means that reconstructing the energy of GCR primaries out of this method is not possible, 
since showers of energies that differ even by a factor 10 can give rise, in some cases, to 
the same position of X^ax- An example from our simulations of the sensitivity of X^ax 
fluctuations to the GCR primary composition is shown in the left panel of fig. [31 Not 
only the position of maxima but also the shape and broadness of the fluctuation profiles 
change according to the GCR masses, as can be seen by comparing the cases of p and Fe 
induced EAS at a fixed energy (10^^ eV for the cases considered in the figure). 

Besides the vertical spreads of particle distribution, also lateral spreads can be cal- 
culated. In particular, in case of EM particles, the lateral spread is due to the effect of 
multiple Coulomb scattering, while, in case of ps, the lateral spread is mainly correlated 
with the pt distribution of charged tt and K, from the decay of whom ps originate (we 
are not mentioning here ps from the decay of short-lived charmed particles, which are 
also considered in FLUKA and other codes). Detailed maps of e and p fluences in the 
atmosphere, down to the Earth's surface, which allow to appreciate both the average 
vertical profile and the lateral spread of each component, as calculated by FLUKA for Fe 
induced showers, can be found in [TO] . 
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